**Figure 5: Counterfactuals**

clear all
use Clean_data, clear

gen femalexcomp= female * inv_competNOQUOTA


** Observed Gender Gap**
gen raw_gg = .

forvalues j = 1983/2019 {
    reg z_total_points i.female##i.year z_GPA if degree == "M" & test == 1 & GPA != . & foreign == 0, robust
    lincom _b[1.female] + _b[`j'.year#1.female]
    
    replace raw_gg = r(estimate) if year == `j'
}

egen mean_merits_points = mean(merits_points), by(year degree) 
egen sd_merits_points = sd(merits_points), by(year degree) 
gen z_merits_points = (merits_points - mean_merits_points) / sd_merits_points


** Gap in Merits Points***

gen GPA_gg = .

forvalues j = 1983/2019 {
    reg z_merits_points i.female##i.year if degree == "M" & test == 1 & GPA != . & foreign == 0, robust
    lincom _b[1.female] + _b[`j'.year#1.female]
    
    replace GPA_gg = r(estimate) if year == `j'
}



**Low Competition***

gen low_comp_gg = .

local min_val = 1e10  

quietly reg z_total_points female femalexcomp z_GPA i.year ///
    if degree == "M" & test == 1 & GPA != . & foreign == 0, robust  

forvalues j = 1983/2019 {
    if inlist(`j', 1987, 1988) continue

    sum inv_competNOQUOTA if year == `j' & degree == "M" & z_GPA != . & foreign == 0, meanonly
    local Pval`j' = r(mean)

    lincom _b[female] + _b[femalexcomp] * `Pval`j''
    local effect = r(estimate)

    // Ensure `effect` is valid before comparing
    if !missing(`effect') & `effect' < `min_val' {
        local min_val = `effect'
    }
}

if `min_val' < 1e10 {
    replace low_comp_gg = `min_val'
}

**High Competition***

gen high_comp_gg = .

local max_val = -1e10  

quietly reg z_total_points female femalexcomp z_GPA i.year ///
    if degree == "M" & test == 1 & GPA != . & foreign == 0, robust  

forvalues j = 1983/2019 {
    if inlist(`j', 1987, 1988) continue

    sum inv_competNOQUOTA if year == `j' & degree == "M" & z_GPA != . & foreign == 0, meanonly
    local Pval`j' = r(mean)

    lincom _b[female] + _b[femalexcomp] * `Pval`j''
    local effect = r(estimate)

    // Ensure `effect` is valid before comparing
    if !missing(`effect') & `effect' > `max_val' {
        local max_val = `effect'
    }
}

if `max_val' > -1e10 {
    replace high_comp_gg = `max_val'
}



*Plot*
twoway line raw_gg GPA_gg low_comp_gg high_comp_gg year, ///
    legend(pos(6) col(4) label(1 "Observed") ///
                          label(2 "Only GPA") ///
                          label(3 "Lowest Competition") ///
                          label(4 "Highest Competition")) ///
    lpattern(solid longdash shortdash dash_dot) /// Define line types
    lcolor(black black black black) 
